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Abstract 



We present a systematic algebraic and numerical investigation of the in- 
' stantaneous Bethe-Salpeter equation. Emphasis is placed on confining inter- 

action kernels of the Lorentz scalar, time component vector, and full vector 
^ \ types. We explore stability of the solutions and Regge behavior for each of 

^ these interactions, and conclude that only time component vector confinement 

leads to normal Regge structure and stable solutions. 



1 Introduction 



The Bethe-Salpeter equation follows from the general principles of quantum 
field theory [Q. The instantaneous form of the Bethe-Salpeter equation, known as 
the Salpeter equation [Q, avoids difficulties related to the relative time degree of 
freedom, and is believed to provide a firm framework for the discussion of bound 
state problems. Almost all knowledge of the solutions to the Salpeter equation is 
restricted to that sector where at least one of the constituent masses is large, in which 
case use of the so called reduced Salpeter equation is justified. Relatively little work 
has been done on the algebraic properties 0, ||, Q of the full Salpeter equation, and 
even less on its numerical solution. A few years ago Lagae proposed a formalism 
which shows considerable promise for the systematic investigation of the full Salpeter 
equation. He also examined several confinement models [§], and concluded that the 
confining potential is not a scalar on the basis of non-linear Regge behavior of the 
equal mass solutions. More recently, Miinz et. al. investigated a linear Lorentz 
scalar or alternatively a time component vector confining kernel combined with an 



effective interaction proposed by 't Hooft from instanton effects in QCD |[T0|. They 
found that there was no convincing parametrization for the confining kernel that 
leads to linear Regge trajectories and also yields spin orbit terms of the correct 
sign. They have also concluded that a scalar confining kernel does not lead to 
stationary solutions for higher angular momenta or small constituent masses. This 
conclusion inspired Parramore and Piekarewicz to perform a stability analysis of 
the variational solutions to the Salpeter equation in the pseudoscalar channel [plj 



They concluded that time component vector confinement is stable with respect to 
the increase in number of basis states, but they found the existence of imaginary 
eigenvalues for scalar confinement. One may argue though [|12] that these authors 
have not considered the norm of the solutions in their analysis. 

In this paper we extend Lagae's method, correcting a small but important al- 
gebraic error and exploring for the first time the full vector interaction kernel. We 
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consider the nature of the solutions to the full Salpeter equation for interaction 
kernels of the time component Lorentz vector, scalar, and full vector types. The 
variational (Galerkin) method is used to investigate the reality of eigenvalues and 
stability of the solutions with the above Lorentz interactions. We have taken the 
norm of the solutions into account in our analysis, since only states with positive 
norm have direct physical significance. We have also investigated the Regge be- 
havior for the above mentioned kernels, and extended the analysis done in |^ to 
heavy-light systems. We also find that only the time component vector interaction 
leads to stable variational solutions and has normal Regge behavior. 

This paper is organized as follows. In Section ^ we review some general properties 
of the Salpeter equation and the reduction to a system of coupled radial equations. 
The numerical results are discussed for each of the Lorentz kernels in Section |^, and 
our conclusions are summarized in Section ^. Appendix ^ contains the complete 
coupled radial equations for the kernels considered in this paper, as well as discussion 
of important limiting cases. The techniques of our numerical solution are covered 
in more detail in Appendix |B|. 

2 Reduction to radial equations 

In this section we briefly review Lagae's elegant formalism ^ for the reduction of the 
full Salpeter equation to a system of equations involving only radial wave functions. 
Our purpose here is to establish notation, and to correct a small inconsistency in the 
Lagae's derivation which, for example, leads to a couple of sign errors in the radial 
equations (5.7) of J^. The correction does not affect the equal mass case, which was 
the subject of Lagae's numerical examples 

We start from the Salpeter equation for a fermion-antifermion system in the CM 
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frame of the bound state, 

Mx = H,x -XH2 + J (^^(k - k')(A^rix'r2A'_ - A^rixTsA^) , (1) 

with notation / = /(k), /' = /(k'). In the above equation V{k — k') is a scalar 
function with Fourier transform V{r) in the case of a Lorentz vector kernel, and 
—V{r) in the case of a Lorentz scalar kernel. The Salpeter amplitude x describing 
a mesonic bound state \B) is defined by 

X(k) = (01^^1(1^)4(1^)15) , (2) 

Fi and F2 are 4x4 matrices corresponding to an interaction kernel with Lorentz 
structure 7°Fi ® 7°F2, while the i^j's and A^'s are generalized Dirac Hamiltonians 
and energy projection operators, given by 

H,{k) = A(k)a-k + S,(k)/5 , (3) 
^ 2E,{k) ' 



with Ei(k) = -y/y4j(k)2 + Bi{ky. In this paper we restrict ourselves to constituent 
quarks of masses rrii, so that 

A(k) = k, (5) 
B,{k) = , (6) 



E,{k) = ^m^ + k^. (7) 

Using properties of projection operators it can be easily shown that the Salpeter 
amplitude satisfies the constraint condition 

|ix + X^ = 0. (8) 



Taking this into account, the norm of the Salpeter amplitude ^ can be written 
as 

r ^3^/ r ff, 1 

(9) 



Ixll = TTT^Tr 



(2vr) 
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and is related to the normalization of bound states as 

IWI^=pi)3(B|B>. (10) 
Using (|I|) inside of d) one obtains 

This equation will be used for obtaining radial equations from the variational prin- 
ciple as outlined in . 

Now, if we expand the Salpeter amplitude as 

X = £o + (^ipi + Afo ■ cr + Afi ■ Pi(T , (12) 

using 16 Hermitian matrices whose squares are unity (1, pi, cr, picr) as defined in 
it is then easily seen that the constraint (§) can be satisfied by expressing the 16 
components of x (-^'s and AT's) in terms of eight functions (Li, L2, Ni, N2) in the 
following way: 

£0 = ^^^(k-Ns) , 
(^i = S^Li , 
C2 = iCgL2 , 
£3 = -C<^(k-Ni) , 

A^o = 5,L2k + 2C^(kx N2) , (13) 
Afi = S^k{k ■ Ni) - Cgk X (k X Ni) , 
Af2 = ^[C,k(k-N2)-5<^kx (kxN2)] , 
Afs = -QLik-25,(kx Ni) . 



Here we have used notation 



= sin (f) , = cos (f) , (14) 
Sg = sin 6* , Cg = cos 9 , (15) 
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with angles and 6 defined as 



+ 



while the 0i's are defined through 



cos ( 



A 



sm ( 



(16) 



(17) 



At this point we have departed from Lagae in one small detail. Namely, we have 
redefined the function Li from equation (4.10) in 0, so that L'j'"'' = ^]^^°-9ae ^ Th.e 
reason for doing so is that now, using ([T^) and (|13D in the expression for the norm 
(^), one can obtain 
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(2^) 



:[L;(k)Li(k) + Lt(k)L2(k) + N;(k)-Ni(k) + Nt(k)-N2(k)] , (18) 



which is equation (4.13) from [^]. Using Lagae's definition of Li would lead to minus 
signs in front of terms L2(k)Li(k) and L^(k)L2(k) in (p^). This small inconsistency 
of equations (4.10) and (4.13) from |^ leads to some incorrect signs in the final form 
of the radial equations for states with parity P = (—1)"'+^. In the equal mass case 
the terms with incorrect signs vanish, so that the numerical results obtained in 
are not affected. 

Now we proceed to obtain the radial equations. We first express Li and 
in terms of spherical harmonics and vector spherical harmonics (for an extensive 
discussion of generic wave functions and the identification of the quantum numbers 
of the bound states the reader is again referred to 0), so that 

L,(k) = L,{k)YjM{k) , (19) 
N,(k) = Ar,_(A:)Y_(k) + Ar,o(A;)Yo(k) + A^,+ (A;)Y+(k) , (20) 

where Y_, Yq, and Y+, stand for Yjj„ia/, Yjjm, and Yjj+im, respectively. We 
also introduce the functions rii^ and ni_, defined as 



(21) 



ni+ 














— V jj 
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with 

Using these definitions, together with properties of spherical and vector spherical 
harmonics, we find from fllSD 



Co 


= Sen2-YjM , 








= S^LiYjM , 






C2 


= 'i-CgL2YjM , 






C3 


= —C(j,niJYjM , 






ATo 


= C^n2+Yo + {uC^N2o + I^SeL2)Y^ + 


if^C^N2o - 


pSeL2)Y+ , 




= CeNioYo + ii^iS^rii. + uCeni+)Y_ ^ 


- {-vS^ni_ 


+ ^LCeni+)Y. 




= i[S^N2oYo + {vS^n2+ + ^iCen2-)Y^ 


+ {^S^n2+ 


- vCen2-)Y^ 




= -Seui+Yo + (-^C^Li - uSgNio)Y. 


_ + {vC^Li 


- ^lSeN^o)Y, 



(24) 



In the above formulas everything is expressed in terms of radial functions (e.g. 
Li = Li{k), and so on). Let us briefiy review the quantum numbers of the states 
that these radial functions represent (parity P, charge conjugation C, and ^'^"*"^Lj): 

Li,L2 P= (-1)^+1 C={-iy 'Jj 

iVio,iV2o P= (-1)^+1 C={-iy+' 'Jj (25) 

ni+,n2+,ni.,n2_ p = {-iy c = {-iy \j±i).j 

Substituting (|2^) in the expression for the norm (0), and using the angular 
integrals summarized in [0, we find 

WxW = 4/ [lil2 + l;l^ + n*,N2o + n;,n,o 

Jo [Zn)-^ 

+ nlj^n2+ + n*2^ni+ + n{__n2- + n*2^niJ\ . (26) 
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Similarly, for the kinetic energy part of ( |Tl| ) we get 



[L^Li + L;L2 + A^ioA^io + iV2o^2o + nl^ni+ + n;+n2+ + nl_ni_ + rf^.n^.] . (27) 



Each type of kernel must be treated separately. Here, for example, we give kernel 
part of (|lT]) for the interaction of the form 7*^ ® 7° (Fi = r2 = 1): 



r r d?\i! , ^ , , \ k'^dk k'^dk' r 

Ll[S^VjS'^L\ + C^(^Vj_i + //Vj+i)C;l; + ^xvC^iVj.^ - Vj^^)S'sN[,] 
+ L;[CeVjC'eL', + Se{fi'Vj^i + z/Vj+i)5^L^ + fiiySe{Vj.i - Vj+i)C'^N^,] 
+ N*,[CeVjC'sN[, + Seiu'Vj^i + fi^Vj+i)S'sN[, + fiuSe{Vj^i - Vj+i)C'^L[] 
+ N;,[S^VjS'^N^, + C<^(/vVj_i + AiVj+i)C;iV^o + f^^C^Vj.i - Vj+,)S'sL',] (28) 
+ nl4SgVjS'sn[^ + Ce{u^Vj^, + i?Vj+r)C,n\^ + iivCe{Vj^^ - Vj+^)S'^n[J^ 
+ nl^[C^VjC'^n'2^ + 5^(z/Vj_i + /iVj+O^^n^^ + iivS^{Vj-i - K;+i)C;n^_] 
+ nl_[C^VjC'^n\_ + 5^(/iVj_i + //Vj+i)5;n;„ + ^iuS^{Vj-i - Vj+i)C'sn\^] 
+ nl_[SeVjS'sn'^_ + Ce(//Vj_i + v^Vj+r)C,n'^_ + ^z/C,(K7_i - Kz+O^^n;,^] } , 

with Vl defined as 

POO 

VL{k,k') = 8n r^drV{r)jL{kr)jL{k'r) . (29) 
Jo 



At this point we can obtain the radial equations by taking variations of (|26|), 
(p7D and (p8D, with respect to Ll{k), L2{k) . . . .n^^ik). The resulting equations for 
the 7° ® 7° kernel, as well as for the 1 1 and 7^ ® 7/^ kernels, are summarized in 
Appendix |^. 

Of course, one can obtain these equations also by straightforward substitution 
of ([T2[) and (|^) into the Salpeter equation (|l|), and then taking the trace after 
multiplication of the resulting equation with the appropriate matrices. The angular 
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integrals which one needs can be worked out easily using the definition of vector 
spherical harmonics, the expression for yj^^(f2)l^™^(f2), and the general properties 
of the Clebsch-Gordan coefficients, as given in the Appendix C of However, 
Lagae's method reviewed in this section is much more simple and elegant. 



3 Numerical results 

As outlined in Appendix |B|, one solves the radial equations by expanding the wave 
functions in terms of a complete set of basis states, which depend on a variational 
parameter (3. This expansion is then truncated to a finite number of basis states. 
In this way, a set of coupled radial equations can be transformed into a matrix 
equation, 7^-0 = Mi/j. The eigenvalues M of the matrix 7i will depend on (3, and by 
looking for the extrema of M{l3), one can find the bound states. If the calculation 
is stable, increasing the number of basis states used will decrease the dependence 
of the eigenvalues on (3. The regions of P with the same eigenvalues should thus 
enlarge. 

A stability analysis of the variational solutions for the pseudoscalar states has 



been recently performed in [|ri[] . Using the fact that the Salpeter equation can be cast 
in a form identical in structure to a random-phase- approximation (RPA) equation, 
the authors of have employed the same formalism developed by Thouless in 



his study of nuclear collective excitations ||Tj] , to perform a stability analysis of the 



Salpeter equation with Lorentz time component vector and scalar confining kernels. 
They find the presence of instability, manifested by the appearance of imaginary 
eigenvalues, in the case of scalar confinement. On the other hand, they find no such 
evidence in the case of time component vector confinement. 

Since matrix Ti is not symmetric, its eigenvalues are not guaranteed to be real. 
However, as noted in 0, the reality of eigenvalues follows from the reality of 
the norm and of the right hand side of (0) when Fi and r2 are hermitian, unless 
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the norm is zero. Physically acceptable solutions must have positive (and non-zero) 
norm. One may argue [0 that the authors of have not taken this into account. 
Therefore, we find it worthwhile to examine the stability of the variational solutions 
to the Salpeter equation by taking the norm into account, i.e. by rejecting states 
with negative or zero norm. 

In the following we consider the above issues for three different kernels, 7° (S> 7°, 
7^ ® IpLi and 1 ® 1. 

3.1 7^ kernel [time component Lorentz vector] 

The case of the 7° ® 7" confining kernel {y{r) = ar) was found |11] to be stable 



with respect to an increase in the number of the basis states, even with very small 
quark masses. In order to verify this, we have performed a similar calculation for the 
pseudoscalar case using a = 0.2 GeV^ and zero mass quarks, with as 

many as 50 basis states. Results for the lowest three physical states (with positive 
norm) are shown in figure |l|. One can see that, as number of basis states increases, 
plateaus with the same eigenvalues enlarge and there is no sign of instability. The 
same calculation was performed with states of higher angular momentum (even 
as high as J = 20), and again results are the same. Therefore, we confirm the 



conclusion reached by the authors of |11| that time component vector confinement 



is well behaved and suitable for a variational solution. As an additional check of 
our programs we have reobtained all the numerical values for the time component 
vector confinement with several different quark masses, which are given in Table 1 
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In Figures || and ^ we show the leading and the first few daughter Regge tra- 
jectories in the case of P = — C mesons in the light-light and heavy-light systems. 
As expected, we obtained slopes of for the light-light, and ^ for the heavy-light 
systems. 

In Figure H we plot radial wave functions for the lowest lying S and P waves in 
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coordinate space. We have used V{r) = ar {a = 0.2 GeV"^) and zero mass quarks. 
The wave functions were normahzed so that = 

Finally, it is a well known fact that the Salpeter equation does not reduce to the 



Dirac equation in the limit where one quark mass becomes infinite |T^. What is 
needed is an interaction which allows the existence of single pair terms. These terms 
arise from kernels involving crossed ladder diagrams. Therefore, it is interesting to 
compare the exact solution of the Dirac equation with the time component vector 
Coulomb potential (y{r) = — ^) with the solution of the Salpeter equation in the 
heavy-light limit. In the Figure^ we plot the Coulomb energy for these two equations 
as a function of the light quark mass, and for three different values of k. For example, 
for the light quark mass of nii = 0.3 GeV and for n = 0.75, the effect is of the order 
of magnitude of 10 MeV. We note that for small mi this system becomes spatially 
large and weakly bound. 

3.2 10 1 kernel [Lorentz scalar] 

As far as stability of the scalar confinement is concerned, things are completely dif- 
ferent. The authors of [|TT| claim that imaginary eigenvalues occur as they increase 



number of basis states, even for large values of the quark mass. For the particular 
choice of mi = m2 = 0.9 GeV and a = 0.29 GeV"^ they find that increasing the 
number of basis states from 20 to 25 leads to the first occurrence of imaginary eigen- 
values (for the pseudoscalar state). Using these parameters we find no imaginary 
eigenvalues, even with as many as 50 basis states. However, as shown in Figure ^ 
scalar confinement in the full Salpeter equation does have a stability problem. As 
soon as basis states having a large enough momentum components are included into 
calculation instabilities occur. As seen in Figure |^, if one includes only 25 basis 
states, the three lowest states with positive norm have all well defined plateaus in 
variational parameter (3. But as soon as we go from 25 to 35 basis the third state 
with positive norm develops instability. Of course, for the ground state this insta- 
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bility occurs later. In order to see the magnitude of this effect we have enlarged the 
scale for the behavior of the ground state plateau in the Figure [7[ As one can see, 
increasing the number of of basis states from 39 to 40 leads to a decrease of energy 
in one small region of (3. Naturally, these problems occur much earlier with smaller 
quark masses. To gain some insight into the nature of this instability, we show in 
Figure ^ the behavior of the wave functions as the instability occurs. We have taken 
the same parameters as for the previous two figures, and chosen (3 = 1.538 GeV (in 
the middle of the instability from the Figure As one might expect, with 39 basis 
states the radial wave functions still preserve behavior characteristic of an S'-wave, 
which is lost for the solutions with 40 basis states. It is also interesting to note that 
(3 near the edge of the unstable region (Figure |^) yields wave functions appear to 
be identical to the ones in the middle of it. 

As far as imaginary eigenvalues obtained from the Salpeter equation with scalar 
confinement are concerned, we have indeed found these with very small quark 
masses. However, states with such eigenvalues always have zero norm, and have 



to be rejected. Therefore, in this case we do not agree with |Tl[, and support con- 
clusions reached by 0, |15|] that physical states will have positive norm and positive 
energies. However, we do agree with |ll[ that the Salpeter equation with scalar 
confinement does have a stability problem. Let us also briefly mention that if one 
squares the ground state energy from Figure |^, one obtains 6.750 GeV^, a number 
to which the calculation of ||Tl|] converges before instability occurs (Table 2 in |1TT| ). 
We do not want to speculate on the reasons why calculations done in are much 
less stable than ours (for example, with parameters used for Figures |^ and we have 
not found imaginary eigenvalues even with 50 basis states, while there they occur 
already with 25 basis states). We might point out though that the pseudo-coulombic 
basis functions we used here are much more suitable for the description of hadronic 



systems than harmonic oscillator basis functions used in [11 



We also confirm that there are problems with Regge trajectories with scalar 
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confinement in the Salpeter equation as already found in H, O]. Witli large quark 



masses and small number of basis states we show in Figure |^ that Regge trajectories 
are not linear. The situation is not improved even in the heavy-light limit, as can 
be seen in Figure |1^. We would like to point out that this can be easily understood 
from the so called "no-pair" equation. The Salpeter equation and its reduced version 
in the heavy-light limit are the same as the no-pair equation |T^, and for the no-pair 
equation with scalar confinement it was shown analytically |TP| that linear Regge 



behavior is lost. If the linear Regge behavior is lost in the heavy-light limit, one 
cannot expect that it will be restored when both quarks have finite mass. 

Finally, let us just mention that we have also investigated the stability of a 
mixture of time component vector and scalar confinement, i.e. for the kernels of the 
type 

X7°(g)7° + (l-x)l(g)l . (30) 

This type of confining kernel (with x = 0.5), together with a one gluon exchange 
kernel, was recently used in for the investigation of the weak decays of B and 
D mesons. In order to illustrate this type of confinement we show in the Figure 



n] what happens as x goes from 0.49 to 0.51 (with zero mass quarks and pure 
confining potential). Obviously, in this case solutions are stable only if x > 0.5. 
With the addition of a short range Coulomb potential, variational solutions exist 



also for X = 0.5. A similar conclusion was also reached in 11 



3.3 7^ kernel [full Lorentz vector] 

Full vector confinement behaves even worse than scalar confinement, as far as the 



variational method is concerned. As one can see in Figure |T2|, the calculation is not 
stable even with quark masses as high as 5.0 GeV. Increasing the number of the 
basis states only makes things worse, as well as decreasing the quark mass. It is 
also interesting to note that this kernel exhibits similar problems even with a pure 
Coulomb potential, which was, on the other hand, found to be stable with scalar 
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and time component vector kernels. Imaginary eigenvalues with this type of kernel 
are quite common, but again we emphasize that all such solutions have zero norm 
and must be rejected. 



4 Conclusions 

In this paper we have corrected a small inconsistency in Lagae's derivation of the 
radial equations for the full Salpeter equation 0, and extended his analysis to the 
case of a full Lorentz vector kernel. We have concentrated here on the nature of 
variational solutions to the full Salpeter equation with a linear confining potential 
and three different types of kernels: time component vector, full vector and Lorentz 
scalar. In each case we have examined the stability of variational solutions, and, 
when possible, the Regge structure in the equal mass and the heavy-light cases. Our 
results support previous conclusions that scalar confinement yields unstable varia- 
tional solutions [O] and non-linear Regge trajectories H, O], even in the heavy-light 



limit. On the other hand, the variational solutions for the time component vector 
confinement are stable, and give linear Regge trajectories with the expected slopes 
for both the equal mass and heavy-light cases. In addition we have found that vari- 
ational solutions for full vector confinement are even more unstable than the ones 
for scalar confinement. We emphasize that our analysis took into account the norm 
of the solutions. Our numerical results completely support the theoretical conclu- 



sions of p, |T5| that eigenvalues for the states with physical norms are always real. 
Imaginary eigenvalues do appear, and are quite common for full vector confinement, 
but these solutions always have zero norm. 

The occurrence of non-linear Regge trajectories in scalar confinement is an ex- 



tension of the results of Gara et. al. (for the reduced Salpeter equation) and 



Lagae (for the full Salpeter equation) in the equal mass case. We have found that 
the same problem persists even in the heavy-light case. The origin of this effect can 
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be explained analytically from the so called "no-pair" equation in coordinate space 

Even though the time component vector interaction behaves perfectly in the full 
Salpeter equation, and leads to the expected linear Regge trajectories, it cannot be 
directly applied in any realistic meson model, because it conflicts with QCD. The 
most evident example of this conflict is its wrong sign of the spin orbit splitting. 
It is important to point out that the relativistic flux tube model reduces to a time 
component vector interaction for S'-waves, but nevertheless yields the correct spin- 



orbit splitting relativistic correction pT |. 
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APPENDICES 



A Radial equations 

In this appendix wc give the complete set of radial equations for the kernels consid- 
ered in this paper. These equations represent the general case of a quark with mass 
nil and an anti-quark with mass 1712- One has to keep in mind that for J = four 
wave functions vanish, i.e. we have A^io = A^20 = 0, and ni+ = 712+ = 0. 

In the equal mass case these equations somewhat simplify, since one has Ei = E2, 
<f) = (pi = (f)2, and = 0, so that Sg = and Cg = 1. Also, since charge conjugation 
is a good quantum number in the equal mass case, the four P = (—1)"'+^ state 
equations split into two systems of two equations, one corresponding to C = (— 1)'^ 
(involving Li and L2), and the other corresponding to C = (—1)'^"'"^ (involving A^io 
and A^2o)- 

The heavy-light limit (m2 — > 00) is obtained by setting £^2 — ^ ?7i2, 02 — ^ f , so 
that Se — > Ccf, and Cg — > S^. It is interesting to note that in the heavy-Ught hmit 
(and for kernels 7^^ (8) 7° and 1 (8) 1, but not for 7'* (8) 7^ kernel) physical solutions 
satisfy Li — L2 and A^io = A^2o (for the P — (— l)-^"^^ states), or ni+ = n2+ and 
?7.i_ = 712- (for the P — {—ly states). Therefore, in these cases the system of four 
equations can be reduced to a system of only two radial equations. 

Of course, for any mixture of different kernels, only the kernel parts of the radial 
equations should be added, and the kinetic energy terms are always the same. In 
the 1 1 case, we have introduced an additional minus sign in the kernel, so that 
V{r) has the same form for all three cases considered, e.g. for the Cornell potential 
V{r) = ar - ^. 



16 



A.l 7" (8) 7" kernel 

States with parity P = {-iy+^: 

ML, = [E, + E2]L2+ jTT^iCeVjC'eL', 

Jo [^T^j 

ML, = [E, + E,]L, + J^ ^[S^VjS'^L[ 

MN^o = [E, + E2]N-,o + -^[S^VjS'^N^o 

+ C^{u^Vj_, + /.Vj+i)C;iV^o + i^i^C^{Vj_r - Vj+,)S',L',] 

+ Se{u^Vj_^ + iJ,^Vj+^)S',N[, + iiiySg{Vj_r - Vj+^)C'^L[] 
States with parity P — (—1)"^: 

POO y^dh' 

Mn,+ = [E^ + E2]n2+ + l j-^[C^VjCy,^ 

poo y'^fiy 
POO y^dh' 

Mni_ = [^i + ^2]n2- + | -^[SeVjS'en'2_ 

POO y^rjh' 

Mn2- = [E, + E2]n,.+ [C^VjC'^n',_ 

Jo [^T^) 
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A. 2 1 (g) 1 kernel 

States with parity P = {-iy+^: 

ML, = [E, + E2]L2+ jTT^iCeVjC'eL', 

Jo [^T^r 

- Se{^^VJ^l + i/Vj+i)5^L^ - nuSeiVj^i - Vj+i)C'^N'^^] 
ML, = [E, + E,]L, + J^ ^[S^VjS'^L[ 

MN^o = [E, + E2]N-,o + -^[S^VjS'^N^o 

- C^{u^Vj_, + /xVj+i)C;7V^o - l^i^C^{Vj_, - Vj+,)S'eL'^] 

States with parity P — (—1)"^: 

Mn,+ = [E, + E2]n2+ + l j^l-C^VjC'^n',^ 

POD y'^dh' 

+ Ceiu^Vj^, + ij,''Vj+Mn[^ + nuCe(Vj-i - Vj+i)S'^n[_] 

roD y^dh' 

Mni_ = [^i + ^2]n2- + | -j^i-SeVjS'Qn'^^ 

+ C,(//Vj_i + v''Vj+,)C'en'2_ + tJivCe{Vj^i - Fj+O^^^+J 

(•oo y'^dh' 

Mn2- = [E, + E2]n,.+ [-C^VjC'^n',_ 

Jo (,^7rj 
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A. 3 7^ kernel 

States with parity P = (— l)*^^"^: 

/•oo y'^dh' 

ML^ = [E, + E2]L2 + J^ j—^[ACeVjC',L'^ 

+ 2Se{fi^Vj^i + u^Vj+i)S'gL', + 2f,uSe{Vj^i - Vj+i)C'^N^o] , 
ML2 = [E, + E,]L,+J^ ^[^2S^VjS'^L[] , (35) 

MN^o = [E, + E2]N2o + -^[2C^(z/Vj_i+/iVj+i)C;iV^o 

MN20 = [E, + E2]N,, + __[2CeVjC'sN[^] . 
States with parity P = (— 1)"': 

roo y^dk' 

Mn,+ = [E^ + E2]n2+ + -^[2C^VjC'^n',^] , 

Mn2+ = [E, + E2]n,+ + -^[2C,(//Vj„i + //Vj+i)CXi+ 

+ 2fiuCeiVj-i - Vj+i)S'^n[_] , (36) 

roo y'^dk' 

Mm_ = [E, + E2]n2-+ -—^[-2SeVjS'en'2_] , 

Jo (,27rj^ 

roo y^dk' 

Mn2- = [E, + E2]n,_+ [AC^VjC'^n\_ 

Jo {271-)^ 

+ 25^(/iVj_i + z/Vj+i)^;r2;_ + 2iivS^{Vj^^ - Vj+^)C'en\^] . 



B Numerical solution of the radial equations 

As briefly mentioned in Section 0, the easiest way to solve the radial equations 
given in Appendix ^ is to expand the radial wave functions f{k) in terms of some 
complete set of basis functions {enik)}, and truncate the expansion to the first 
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basis states, i.e. 



Af-l 



(37) 



i=0 



In this way a set of n coupled radial equations becomes the nN x nN matrix equation 

Htlj = Mij , (38) 
where i(j is an nN dimensional vector, 



(39) 



\fn J 

The matrix Ti. and its eigenvalues depend on the variational parameter /? character- 
izing the basis functions. However, if the calculation is stable, dependence of the 
solution on /? should reduce as increases. This is manifested by the development 
of plateaus in (3 having the same eigenvalues. 

A basis set that was shown to be very successful in calculations of this sort is 
given by [P] 



e.,(Afc) = (-^)-iV.,/3i(.:A.)^-^fc-#"-^"^^^^^^ 



(40) 



where P^"''^\x) are Jacobi polynomials and 



2T( 



iL 



i\(i + 2L + 2) 



TT 



r(z + L + |) 

The Fourier transform of (^OD is known analytically, 

e,L{P,r) = Ar,i/?i(2/3r)^e-^^Lf^+')(2/3r) 
with L'f\x) being the generalized Laguerre polynomials and 



iL 



(i + 2L + 2)!_ 



(41) 



(42) 



(43) 
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The kinetic energy terms of the Ti matrix (dropping the dependence of the basis 
states on 

roo 

/ k^dke*L{k)Vk^ + m^ejLik) , (44) 

<^ 

can be efficiently calculated using Gauss- Jacobi quadrature formula after performing 
a change of integration variable to x = ^2^^- Keeping in mind the definition (|29|), 
the kernel part of 7i will in general include terms like 

O roo roo roo 

- / k^dk / k'^dk' / r^dre*L{k)F{k)]L'{kr)V{r)jL>{k'r)G{k')e,L{k') , (45) 

TT JO Jo Jo 

where F and G are some functions of k and k', respectively. It is quite difficult to 
approximate these integrals by applying standard numerical quadrature methods, 
since the range of integration is infinite and spherical Bessel functions ji{kr) are 
rapidly oscillating. Besides, we would have to perform the quadrature for each 
value of r where V{r)jLr{kr)jLi{k'r) is needed. Therefore, it is much better to 
expand F{k)eiL{k) in terms of basis functions {ei'L'} i.e. 

N'-l 

F{k)e,L{k) ^ cS?^'^^e..i,(^) , (46) 

i'=0 

where A^' is the number of basis states of {eiiy} used (similarly we have to expand 
G{k')ejL{k')). Now, the radial Fourier transforms, 

^7]- roo 

e^Lir) = t^——T / k^dkjLikr)e^Lik) , (47) 
(27r)2 Jo 

allow us to perform integrations over k and k' analytically, so that (^) becomes 

E'E'4!"''"v.if45^'", (48) 

i'=0 j'=0 

where, 

4f^'^) = r edke*,L,{k)F{k)e,L{k) , (49) 
Jo 



V^^P = / r'dre:,^,{r)V{r)ej,L'{r) . (50) 



oo 
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Again, the integrals involved in the calculation of c\n' ' can be efficiently eval- 
uated using the Gauss- Jacobi quadrature formula, while matrix elements of V{r) 
were calculated analytically [|Ty[ for the short range Coulomb and linear confining 
potential that were used in this paper. 

In practice, we choose L in such a way that L = J for the states with parity 
P = (—1)"^+^, and L = J — 1 for the states with parity P = (— 1)"^ (unless J = 0, 
when we take L = 1. From the radial equations given in Appendix ^ one can 
see that L' can take values J — 1, J, or J + 1, and we have found that generally 
taking A^' = 2N is more than enough to accurately describe all functions of the form 
F{k)e,L{k). 

In order to make sure that numerical calculations were done correctly, we have 
written two independent programs, each of them treating the general (unequal mass) 
case, and also two limiting cases (equal mass and heavy-light limit) separately, and 
for all kernels considered in the text. Besides consistency of the two programs of 
ours, and consistency of the two special cases with the general case, we have also 
checked our results against those of 0, where only equal mass case was considered, 
and of |TT|, where stability of pseudoscalar case was considered. 
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FIGURES 



Figure 1: The lowest three states with positive norm for time component vector 
confinement {V{r) = ar, a = 0.2 GeV"^) in the pseudoscalar case, with mi = m2 = 0. 
Calculations were done using 5, 15, 25 and 50 basis states. 

Figure 2: Regge trajectories for a time component vector confinement {V{r) = ar, 
a = 0.2 GeV^) in the light-light case with equal mass quarks P = (—1)'^"'"^ and 
C = (—1)"'. We have taken nii = m2 = 0, and used 25 basis states. 

Figure 3: Regge trajectories for a time component vector confinement {V{r) = ar, 
a = 0.2 GeV"^) in the heavy-light case with P = (—1)"'+^. We have taken mi = 0, 
7712 = 5.0 GeV, and used 25 basis states in the equations corresponding to the 
heavy-light limit. 

Figure 4: Radial wave functions in coordinate space for the J^*" = 0""*" (^5*0 state, 
Li is the lower, and L2 the upper full line), and J'^*^ = O^"*" (^Pq state, rii^ is the 
upper, and 77-2- the lower dashed line), with a time component vector kernel and 
V{r) = ar {a — 0.2 GeV"^). We have chosen quark masses mi — 1712 — 0. The 
calculation was done with 15 basis states. 

Figure 5: Comparison of the Coulombic ground state energy as a function of light 
quark mass mi, obtained from Dirac the (full lines) and heavy-light Salpctcr equa- 
tion (dashed lines), with a time component vector Coulomb potential {V{r) = — ^). 

Figure 6: The lowest three states with positive norm for scalar confinement (y{r) — 
ar, a — 0.29 GeV"^) in the pseudoscalar case, with mi = m2 = 0.9 GeV. Calculations 
were done using 5, 25 and 35 basis states. 



25 



Figure 7: The lowest pseudoscalar state with positive norm for scalar confinement 
(V{r) = ar, a = 0.29 GeV^) with mi = ms = 0.9 GeV. Calculations with 39 
(dashed line) and 40 (full line) basis states are shown. This illustrates the onset of 
instability for this example of scalar confinement. 

Figure 8: Pseudoscalar ground state radial wave functions Li and L2 in coordinate 
space with a scalar confining kernel and V{r) = ar {a = 0.29 GeV"^). Quark masses 
were mi = m2 = 0.9 GeV. Calculations with 39 (dashed line) and 40 (full line) 
basis states are shown, and variational parameter was (3 = 1.538 GeV. 

Figure 9: Regge trajectories for scalar confinement {V{r) = ar) with equal mass 
quarks and P = (— 1)"'"''^ and G = (— 1)"'. We have taken mi = 1112 = 1 GeV, 
a = 0.2 GeV^, and used 25 basis states. 

Figure 10: Regge trajectories for scalar confinement {V{r) = ar, a = 0.2 GeV"^) 
in the heavy-light case with P = (—1)"^+^. We have taken mi = 0.1 GeV, m2 = 
5.0 GeV, and used 25 basis states in the equations corresponding to the heavy-light 
limit. 

Figure 11: The lowest state with positive norm for the mixture of time component 
vector and scalar confinement as given in ( PD| ) {V{r) = ar, a = 0.2GeV^, mi = 
m2 = 0). Calculations with x = 0.51 (full line), x = 0.50 (dashed line), and 
X = 0.49 (dotted line), correspond to 51%, 50%, and 49% of the time component 
vector kernel. We have used 15 basis states. 

Figure 12: The lowest state with positive norm for full Lorentz vector confinement 
iV{r) = ar, a = 0.2 GeV^) in the pseudoscalar equal mass case, with the choice 
of mi = m2 = 5.0 GeV. Calculations with 5 (dashed line) and 15 (full line) basis 
states are shown. 
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